# prepared by Sam. feel free to consult (sirmaxford@gmail.com).
import matplotlib.pyplot as plt; import pandas as pd; import numpy as np
from scipy.optimize import curve_fit
from scipy.interpolate import interp1d
from scipy.signal import savgol_filter
from scipy.interpolate import UnivariateSpline
%matplotlib inline
font = {'family': 'serif',
'color': 'darkred',
'weight': 'normal',
'size': 14,
}
columns = ['r', 'b', 'atp', 'md', 'v_av', 'v_sd']
df = pd.read_csv('vel_changes.csv', names=columns)
df[0:10]
r_change = df[0:10]
r_change['v_av']
x = r_change['r']
y = r_change['v_av']
xx = np.linspace(x.min(), x.max(), 1000)
# interpolate + smooth
itp = interp1d(x,y, kind='linear')
window_size, poly_order = 101, 3
yy_sg = savgol_filter(itp(xx), window_size, poly_order)
# or fit to a global function
def func(x, A, B, x0, sigma):
return A+B*np.tanh((x-x0)/sigma)
fit, _ = curve_fit(func, x, y) #, bounds=(0.0,[1.0,0.3,0.6]))
yy_fit = func(xx, *fit)
#s = UnivariateSpline(x, y, s=2)
#xs = np.linspace(x.min(), x.max(), 100)
#ys = s(xs)
$f(x) = A + B \tanh\left(\frac{x-x_0}{\sigma}\right)$
#plt.figure(figsize=(10,8))
plt.figure(figsize=(14,10), dpi=500)
#plt.plot(xx, yy_sg, 'k', label= "Smoothed curve")
#plt.plot(xs, ys)
plt.plot(xx, yy_fit, 'c--', linewidth=3.5, label= "Curve fit")
plt.errorbar(x,y,yerr=r_change['v_sd'], ecolor='b', capsize=7, \
mec='blue', color='green', marker='D', linewidth=2, \
markersize=7, label='Mean Velocity Deviation: ATP = 2000, MD = 3000')
plt.xlabel('Specie Ratio', fontdict=font); plt.ylabel('Mean Velocity', fontdict=font)
plt.title('Actin Filament Movement'); plt.legend(loc='upper left'); plt.grid()
df[10:25]
atp_change = df[10:25]
x = atp_change['atp']
y = atp_change['v_av']
xx = np.linspace(x.min(), x.max(), 1000)
# interpolate + smooth
itp = interp1d(x,y, kind='linear')
window_size, poly_order = 101, 3
yy_sg = savgol_filter(itp(xx), window_size, poly_order)
# or fit to a global function
#def func(x, A, B, x0, sigma):
# return A+B*np.tanh((x-x0)/sigma)
def func(x, a, b, c):
return a * np.exp(-b * x) + c
fit, _ = curve_fit(func, x, y)
yy_fit = func(xx, *fit)
plt.figure(figsize=(14,10), dpi=500)
#plt.plot(xx, yy_fit, 'c--', linewidth=3.5, label= "Curve fit")
plt.errorbar(x,y,yerr=atp_change['v_sd'], ecolor='b', capsize=7, \
mec='blue', color='green', marker='D', linestyle='-', linewidth=2, \
markersize=7, label='Mean Velocity Deviation: Ratio = 90%, MD = 3000')
plt.xlabel('ATP', fontdict=font); plt.ylabel('Mean Velocity', fontdict=font)
plt.title('Actin Filament Movement'); plt.legend(loc='upper left'); plt.grid()
df[25:30]
md_change = df[25:30]
plt.figure(figsize=(14,10), dpi=500)
plt.errorbar(md_change['md'],md_change['v_av'],yerr=md_change['v_sd'], ecolor='b', capsize=7, \
mec='blue', color='green', marker='D', linestyle='-', linewidth=2, \
markersize=7, label='Mean Velocity Deviation: Ratio = 90%, ATP = 2000')
plt.xlabel('Motor Density', fontdict=font); plt.ylabel('Mean Velocity', fontdict=font)
plt.title('Actin Filament Movement'); plt.legend(loc='upper left'); plt.grid()
conf_names = ['beads', 'x', 'y', 'z']
df_conf = pd.read_csv('Conformation_A001.txt', names=conf_names, delim_whitespace=True)
df_conf_good = pd.read_csv('Conformation_A001_good.txt', names=conf_names, delim_whitespace=True)
df_conf = df_conf.drop(['beads', 'x', 'y'], axis=1)
df_conf_good = df_conf_good.drop(['beads', 'x', 'y'], axis=1)
df_conf = df_conf[13:]
df_conf_good = df_conf_good[13:]
df_z = df_conf.iloc[0::13, :]
df_z_good = df_conf_good.iloc[0::13, :]
df_z.shape
tip_names = ['time', 'x_tip', 'y_tip']
df_tip = pd.read_csv('TipXY_A001.txt', names=tip_names, delim_whitespace=True)
df_tip_good = pd.read_csv('TipXY_A001_good.txt', names=tip_names, delim_whitespace=True)
df_time = df_tip.drop(columns=['x_tip','y_tip'])
df_time_good = df_tip.drop(columns=['x_tip','y_tip'])
df_time = df_time[1:]
df_time_good = df_time[0:]
df_time_good.shape
Dx_tip = np.diff(df_tip['x_tip']); Dy_tip = np.diff(df_tip['y_tip'])
Dx_tip_good = np.diff(df_tip_good['x_tip']); Dy_tip_good = np.diff(df_tip_good['y_tip'])
dist=np.sqrt((Dx_tip**2)+(Dy_tip**2))
dist_good=np.sqrt((Dx_tip_good**2)+(Dy_tip_good**2))
tym=0.01
vel_tip = dist/tym
vel_tip_good = dist_good/tym
df_vel_tip = pd.DataFrame(data=vel_tip)
df_vel_tip_good = pd.DataFrame(data=vel_tip_good)
df_vel_tip.shape

fig, ax1 = plt.subplots(dpi=500)
ax1.plot(df_time_good,df_z_good, 'b', label='Actin Z-Variations')
ax1.set_xlabel('Time (s)',fontdict=font)
ax1.set_ylabel('Z-Distance', color='b',fontdict=font)
ax1.tick_params('y', colors='b')
plt.legend(loc='upper left')
ax2 = ax1.twinx()
ax2.plot(df_time_good,df_vel_tip_good, 'r', label='Tip Velocity')
fig.set_size_inches([15, 7])
fig.tight_layout()
ax2.set_ylabel('Tip-Velocity', color='r',fontdict=font)
ax2.tick_params('y', colors='r')
plt.title('Actin Filament Movement',fontdict=font); plt.legend(loc='upper right')
plt.grid(); plt.show()
#df_z_good.to_csv('df_z_good', encoding='utf-8', index=False)
#df_vel_tip_good.to_csv('df_vel_tip_good', encoding='utf-8', index=False)
c_names = ['Change in Z', 'Tip Velocity']
z_and_tip_vel_good = pd.read_csv('z_and_tip_vel_good.csv', names=c_names)
z_and_tip_vel_good.corr()
fig = plt.figure(figsize=(14,10), dpi=500)
ax = fig.gca()
df_z_good.hist(bins=100, ax=ax, color='green')

fig, ax1 = plt.subplots(dpi=500)
ax1.plot(df_time,df_z, 'b', label='Actin Z-Variations')
ax1.set_xlabel('Time (s)',fontdict=font)
ax1.set_ylabel('Z-Distance', color='b',fontdict=font)
ax1.tick_params('y', colors='b')
plt.legend(loc='upper left')
ax2 = ax1.twinx()
ax2.plot(df_time,df_vel_tip, 'r', label='Tip Velocity')
fig.set_size_inches([15, 7])
fig.tight_layout()
ax2.set_ylabel('Tip-Velocity', color='r',fontdict=font)
ax2.tick_params('y', colors='r')
plt.title('Actin Filament Movement',fontdict=font); plt.legend(loc='upper right')
plt.grid(); plt.show()
#df_z.columns = ['Change in Z']
#df_vel_tip.columns = ['Tip Velocity']
#df_z = df_z.reset_index(drop=True)
#df_z.merge(df_vel_tip, how='left', left_on='Change in Z', right_on='Tip Velocity')
#np.array(df_z, df_vel_tip)
#df_z.to_csv('df_z', encoding='utf-8', index=False)
#df_vel_tip.to_csv('df_vel_tip', encoding='utf-8', index=False)
c_names = ['Change in Z', 'Tip Velocity']
z_and_tip_vel = pd.read_csv('z_and_tip_vel.csv', names=c_names)
z_and_tip_vel.corr()
fig = plt.figure(figsize=(14,10), dpi=500)
ax = fig.gca()
df_z.hist(bins=100, ax=ax, color='green')
from mpl_toolkits import mplot3d
conf_names = ['beads', 'x', 'y', 'z']
df_xyz = pd.read_csv('Conformation_A001.txt', names=conf_names, delim_whitespace=True)
df_xyz_good = pd.read_csv('Conformation_A001_good.txt', names=conf_names, delim_whitespace=True)
df_xyz_1 = pd.read_csv('Conformation_A001_1.txt', names=conf_names, delim_whitespace=True)
df_xyz_2 = pd.read_csv('Conformation_A001_2.txt', names=conf_names, delim_whitespace=True)
df_xyz_3 = pd.read_csv('Conformation_A001_3.txt', names=conf_names, delim_whitespace=True)
df_xyz_4 = pd.read_csv('Conformation_A001_4.txt', names=conf_names, delim_whitespace=True)
fig = plt.figure(figsize=(14,10), dpi=500)
ax = plt.axes(projection='3d')
ax.plot3D(df_xyz['x'], df_xyz['y'], df_xyz['z'], 'magenta')
ax.scatter3D(df_xyz['x'], df_xyz['y'], df_xyz['z'], 'green') #, c=df_xyz['z'], cmap='Greens')
#ax.view_init(90, 0)
ax.set_xlabel('X Tip')
ax.set_ylabel('Y Tip')
ax.set_zlabel('Z Tip')
fig = plt.figure(figsize=(14,10), dpi=500)
ax = plt.axes(projection='3d')
ax.plot3D(df_xyz_good['x'], df_xyz_good['y'], df_xyz_good['z'], 'magenta')
ax.scatter3D(df_xyz_good['x'], df_xyz_good['y'], df_xyz_good['z'], 'green')
#ax.view_init(90, 0)
ax.set_xlabel('X Tip')
ax.set_ylabel('Y Tip')
ax.set_zlabel('Z Tip')
#fig = plt.figure(figsize=(14,10), dpi=500)
#ax = plt.axes(projection='3d')
#df_new = df_xyz.unstack().reset_index()
#df_xyz['x']=pd.Categorical(df_xyz['x']); df_xyz['y']=pd.Categorical(df_xyz['y']); df_xyz['z']=pd.Categorical(df_xyz['z'])
#df_xyz['x']=df_xyz['x'].cat.codes; df_xyz['y']=df_xyz['y'].cat.codes; df_xyz['z']=df_xyz['z'].cat.codes
#ax.plot_trisurf(df_xyz['x'], df_xyz['y'], df_xyz['z'], 'magenta')
#ax.scatter3D(df_xyz['x'], df_xyz['y'], df_xyz['z'], 'green')
fig = plt.figure(figsize=(14,10), dpi=500)
ax = plt.axes(projection='3d')
ax.plot3D(df_xyz_1['x'], df_xyz_1['y'], df_xyz_1['z'], 'magenta')
ax.scatter3D(df_xyz_1['x'], df_xyz_1['y'], df_xyz_1['z'], 'green')
#ax.view_init(90, 0)
ax.set_xlabel('X Tip')
ax.set_ylabel('Y Tip')
ax.set_zlabel('Z Tip')
fig = plt.figure(figsize=(14,10), dpi=500)
ax = plt.axes(projection='3d')
ax.plot3D(df_xyz_2['x'], df_xyz_2['y'], df_xyz_2['z'], 'magenta')
ax.scatter3D(df_xyz_2['x'], df_xyz_2['y'], df_xyz_2['z'], 'green')
#ax.view_init(90, 0)
ax.set_xlabel('X Tip')
ax.set_ylabel('Y Tip')
ax.set_zlabel('Z Tip')
fig = plt.figure(figsize=(14,10), dpi=500)
ax = plt.axes(projection='3d')
ax.plot3D(df_xyz_3['x'], df_xyz_3['y'], df_xyz_3['z'], 'magenta')
ax.scatter3D(df_xyz_3['x'], df_xyz_3['y'], df_xyz_3['z'], 'green')
#ax.view_init(30, 270)
ax.set_xlabel('X Tip')
ax.set_ylabel('Y Tip')
ax.set_zlabel('Z Tip')
fig = plt.figure(figsize=(14,10), dpi=500)
ax = plt.axes(projection='3d')
ax.plot3D(df_xyz_4['x'], df_xyz_4['y'], df_xyz_4['z'], 'magenta')
ax.scatter3D(df_xyz_4['x'], df_xyz_4['y'], df_xyz_4['z'], 'green')
# rotate the axes and update
#for angle in range(0, 360):
# ax.view_init(30, angle)
# plt.draw()
# plt.pause(.001)
#ax.view_init(90, 0)
ax.set_xlabel('X Tip')
ax.set_ylabel('Y Tip')
ax.set_zlabel('Z Tip')